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HIGHLIGHTS 


►  We  develop  a  comprehensive  diffusivity  model  for  gas  diffusion  layers  (GDLs)  in  proton  exchange  membrane  fuel  cells  (PEMFCs). 

►  The  analytical  model  agrees  excellently  with  experimental  results  and  numerical  data  available  in  literature. 

►  The  influences  of  microstructures  of  GDLs  are  extensively  explored. 
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In  gas  diffusion  layers  (GDLs)  of  proton  exchange  membrane  fuel  cells  (PEMFCs),  effective  gas  diffusivity 
is  a  key  parameter  to  be  determined  and  engineered.  Existing  theoretical  models  of  effective  diffusivity 
are  limited  to  one-dimensional  (ID)  regular  fiber  arrays.  Numerical  simulations  were  carried  out  to 
simulate  gas  diffusion  through  more  realistic  fibrous  materials  like  GDLs,  in  which  fibers  are  randomly 
distributed  in  a  two-dimensional  (2D)  plane  or  three-dimensional  (3D)  space,  but  they  could  not  fully 
reveal  the  underlying  mechanisms.  In  this  paper,  we  propose  an  analytical  model  to  predict  the  effective 
diffusivities  of  ID,  2D  and  3D  randomly  distributed  fiber  assembles.  The  present  model  is  established  by 
extending  the  model  of  ID  regular  fiber  alignments  to  ID  random  fiber  arrangements  through  Voronoi 
Tessellation  method,  and  using  the  ID  local  diffusivities  to  determine  the  2D  and  3D  diffusivities  based 
on  mixing  rules.  The  predicted  effective  diffusivities  agree  well  with  experimental  results  and  numerical 
data.  With  the  new  model,  the  influences  of  porosity,  fiber  distribution,  and  fiber  orientation  are 
analyzed  in  this  study. 

©  2012  Elsevier  B.V.  All  rights  reserved. 


1.  Introduction 

Proton  exchange  membrane  fuel  cell  (PEMFC)  is  considered  to 
be  one  of  the  leading  candidates  for  the  power  sources  of  mobile, 
stationary,  and  portable  devices  [1  ].  The  gas  diffusion  layer  (GDL)  of 
PEMFCs  is  a  fibrous  porous  material  with  a  layered  structure,  which 
not  only  provides  the  support  of  the  fuel  cell  membrane,  but  also 
allows  the  transport  of  reactant  products.  For  example,  oxygen 
diffuses  through  the  GDL  from  the  gas  channel  (GC)  to  the  catalyst 
layer  (CL),  where  it  is  combined  with  the  protons  and  electrons 
from  the  anode  to  produce  water  [2],  The  produced  water  could 
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condense  and  even  block  the  porous  GDL.  Therefore,  effective 
diffusivities  of  water  vapor,  oxygen,  and  hydrogen  strongly  affect 
the  PEMFCs  performance. 

The  movement  of  gas  molecules  caused  by  concentration 
difference  in  a  porous  medium  is  known  as  diffusion  [3].  It  takes 
place  when  the  concentration  of  the  molecules  is  higher  in  one 
region  than  the  other.  Gas  molecules  will  not  stop  migrating  until 
there  is  an  equalized  concentration  configuration  throughout  the 
carrier.  The  moving  paths  of  molecules  during  diffusive  motion 
process  are  random,  but  the  most  preferred  migration  of  molecules 
will  be  in  the  direction  of  decreasing  concentration.  Generally,  gas 
diffusion  through  a  medium  can  be  phenomenologically  described 
by  Fick’s  law: 

J  =  DeffVC,  (1) 

where,  J  is  the  diffusive  flux,  Def r  is  the  effective  diffusivity  and 
VC  is  the  concentration  gradient.  Fick’s  law  states  that  the  average 
flux  in  the  fibrous  structure  is  directly  proportional  to  the  gas 
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concentration  gradient,  and  the  effective  diffusivity  tensor  De ff 
lumps  all  the  complex  interactions  between  gas  and  fiber.  Accurate 
determination  of  effective  diffusivity,  therefore,  is  essential  to 
characterize  the  diffusion  process. 

The  effective  diffusivity  of  fibrous  materials  can  be  determined 
through  Eq.  (1)  by  experimentally  measuring  the  diffusive  flux  and 
the  directional  concentration  gradient.  Early  in  1940,  through-plane 
diffusion  coefficient  of  a  steel  wool  sample  was  measured  using 
carbon  disulfide  and  acetone  vapor  [4],  In  1984,  Bateman  et  al.  [5] 
employed  NO  gas  to  transfer  through  a  2D  cellulosic  filter  and 
calculated  its  effective  diffusivity.  Recently,  Gibson  et  al.  [6]  applied 
a  dynamic  diffusion  test  cell  method  to  measure  vapor  diffusivity  of 
fibrous  media.  Two  parallel  gas  flows  with  different  water  vapor 
humidity  were  directed  to  a  test  cell,  with  which  vapor  diffused 
through  the  sample  between  the  gas  flows,  and  the  effective  diffu¬ 
sivity  was  determined  by  measuring  the  relative  humidity  of  gas 
flows  leaving  the  cell  [6],  Huang  and  Qian  [7]  modified  the  dynamic 
cell  method  by  providing  a  water  vapor  source  on  one  side  of  the 
sample.  Effective  diffusivities  of  GDLs  were  also  obtained  indirectly 
by  measuring  ionic  conductivity  for  the  soaped  electrolyte  based  on 
the  analogy  between  Ohm’s  law  and  Fick’s  law  [8—10],  Using 
a  dynamic  diffusion  cell,  LaManna  and  Kandlikar  [11]  recently 
investigated  the  effects  of  Microporous  Layer  (MPL)  coatings,  GDL 
thickness,  and  polytetrafluorethylene  (PTFE)  loadings  on  the  effective 
water  vapor  diffusion  coefficient  of  GDLs.  Instead  of  measuring  the 
effective  diffusivity,  many  measured  the  evaporative  moisture  vapor 
resistance,  which  is  inversely  related  to  the  vapor  diffusivity  [12,13]. 

Many  researchers  modeled  the  gas  diffusion  through  fibrous 
materials.  The  simplest  model  assumes  the  fibrous  material  con¬ 
sisting  of  a  bundle  of  tortuous  channels.  So  the  effective  diffusivity 
is  related  to  the  bulk  diffusivity  in  the  void  through  porosity  and 
tortuosity,  given  by  a  normalized  form  [14],  viz: 


where  Db  is  the  gas  diffusivity  in  the  void,  t  is  the  porosity,  and  t  is 
the  tortuosity.  Although  porosity  is  easy  to  be  calculated  or 
measured,  the  applicability  of  Equation  (2)  is  tarnished  by  the 
difficulty  in  accurately  determining  the  value  of  tortuosity  [15].  A 
pore-scale  model  was  presented  to  predict  the  effective  diffusivity 
of  unconsolidated  porous  media  based  on  a  rectangular  represen¬ 
tative  unit  cell,  in  which  the  tortuosity  was  expressed  as  the  ratio  of 
the  diffusive  path  length  to  the  streamwise  displacement  [16]. 
However,  the  regular  geometry  of  the  diffusive  streamlines  in  the 
pore-scale  model  is  over-idealized. 

Apart  from  the  pore-based  models  above,  fiber-based  models 
were  developed.  Shen  and  Springer  [17]  calculated  diffusive 
transport  across  ID  impermeable  cylinders  with  square  packing 
configuration  in  a  rectangular  unit  cell,  and  the  model  of  effective 
diffusivity  was  expressed  as: 

I f  =  1-2V^’  (3) 

which  was  widely  applied  to  evaluate  the  influence  of  water  vapor 
diffusion  on  the  mechanical  properties  of  composites  [18],  Never¬ 
theless,  this  model  did  not  consider  the  varying  width  of  the  gaps 
between  cylindrical  fibers  and  hence  under-estimates  the  effective 
diffusivity.  The  varying  width  of  the  gap  between  cylindrical  fibers 
was  later  considered  by  Li  et  al.  [19]  for  both  square  and  hexagonal 
fibers  arrangement,  but  their  model  over-predicts  the  diffusivity  as 
the  gas  concentration  is  assumed  to  be  constant  at  any  cross-section 
in  the  channel  between  the  fibers.  In  the  PEMFC  literature,  the 
Bruggeman  model  has  attracted  most  of  the  attention  [20],  The 


Bruggeman  model  is  based  on  effective  medium  approximation  and 
is  given  by: 


However,  the  Bruggeman  model  was  derived  for  uniformly 
packed  spherical  particles  rather  than  differently  oriented  cylin¬ 
drical  fibers. 

In  order  to  model  gas  diffusion  through  more  realistic  fibrous 
structures,  a  number  of  researchers  applied  a  variety  of  numerical 
simulation  techniques.  Tomadakis  and  Sotirchos  [14]  performed 
Monte  Carlo  simulations  for  ID,  2D  and  3D  randomly  oriented 
fibers  and  calculated  the  effective  diffusivity.  They  measured  the 
mean  traveling  distances  of  diffusive  molecules  inside  the  numer¬ 
ical  fiber  network,  and  proposed  the  following  model  for  the 
randomly  oriented  fibrous  systems: 


where,  tp  is  the  percolation  threshold  and  is  an  empirical  constant 
determined  by  a  least  squares  fit  to  the  simulation  results.  In 
another  numerical  study,  the  local  effective  diffusivity  of  a  GDL 
medium  was  determined  as  a  function  of  the  local  porosity  and  the 
local  water  saturation  by  a  network  model,  where  the  solid  struc¬ 
ture  was  simulated  as  layers  of  fiber  screens  and  each  layer  was 
shifted  by  a  randomly  selected  in-plane  distance  [21],  In  2008, 
Becker  et  al.  [22]  numerically  reconstructed  a  fibrous  structure 
from  a  3D  tomography  image  of  the  GDL,  and  proposed  the  effec¬ 
tive  diffusivity  as  a  function  of  the  saturation  of  the  GDL.  Later  in 
2011,  Becker  et  al.  [23]  extended  their  work  to  consider  the  effect  of 
Microporous  Layer  (MPL). 

Although  many  analytical  models  have  been  proposed  for 
determining  the  effective  diffusivity  of  fibrous  materials,  they  are 
limited  to  1 D  fiber  arrays.  For  more  realistic  fibrous  materials  where 
fibers  are  2D  or  3D  oriented,  only  numerical  studies  have  been 
conducted,  which  cannot  fully  reveal  the  underlying  mechanisms 
of  gas  diffusion.  In  this  work,  we  propose  a  generalized  analytical 
model  of  effective  diffusivity  of  fibrous  materials.  The  new  model  is 
established  through  first  extending  the  model  of  regular  ID  fiber 
arrays  to  random  ID  fiber  arrays  by  Voronoi  Tessellation  approxi¬ 
mation,  and  then  applying  them  to  2D  and  3D  random  fiber 
assembles  by  mixing  rules. 

2.  Model  generation 

In  this  paper,  fibrous  media  are  assumed  to  be  composed  of 
periodical  unit  cells  representing  the  geometrical  nature  of  the 
microstructures.  The  present  model  is  established  by  extending  the 
model  for  regular  ID  fiber  arrays,  since  although  the  architectures 
of  fibrous  media  vary  from  simple  ID  regular  type  to  complex  2D 
and  3D  cases,  they  can  be  approximately  represented  by  mixtures 
of  ID  fibers  arrays  [24,25],  The  following  assumptions  are  made: 

1.  The  fibrous  matrix  is  made  up  of  straight  and  circular  fibers 
with  relatively  high  porosity. 

2.  All  the  fibers  are  impermeable  and  the  gas  diffusion  takes  place 
only  in  the  void  areas  between  fibers. 

3.  The  fibrous  media  have  relatively  high  porosity  and  the  spacing 
between  fibers  is  much  larger  than  the  mean  free  path  of 
diffusing  species. 

The  simplest  representative  cell  for  ID  fibrous  media  is  regular 
array  of  parallel  fibers  as  shown  in  Fig.  1.  The  diffusivity  in  the  open 
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voids  of  the  cell  is  equal  to  bulk  diffusivity  Db,  which  characterizes 
gas  diffusion  in  bulk  space  without  confinement.  This  representa¬ 
tive  cell  is  composed  of  an  impermeable  fiber  and  the  gas  matrix, 
surrounded  by  the  effective  medium  with  diffusivity  Defr,  which  is 
the  same  as  the  diffusivity  of  the  whole  system.  For  convenience 
and  without  losing  generality,  the  representative  cell  is  assumed  to 
be  a  circle  with  the  same  area  of  the  square  containing  the  fiber.  The 
porosity  c  for  this  arrangement  of,  both  the  unit  cell  and  the  whole 
fibrous  system  can  be  determined  by  the  following  equation: 


Q2  “  RT 


(6) 


where  r  is  the  fiber  radius,  R  is  the  radius  of  the  unit  cell,  a  is  the 
edge  length  of  the  square.  In  a  steady  state,  the  gas  diffusion  across 
the  fiber  array  is  defined  by  the  following  Laplace  equation  [26]: 

V2C(/,0)  =  0,  (7) 


where  1  and  6  are  the  cylindrical  coordinates.  The  appropriate 
solution  of  Eq.  (7)  for  diffusion  normal  to  fibers  reads  [27]: 


!Al  cos  6,  (0  <  /  <  r) 

(B/  +  0cos0,(r</<R),  (8) 

Cl  cos  8,  (R  <  l) 


where  A,  B,  E,  and  G  are  unknown  constants  arising  from  boundary 
conditions. 

The  continuity  of  the  diffusive  gas  at  the  boundary  between  the 
impermeable  fiber  and  the  gas  leads  to: 

Br  +  -  =  Ar  =  0,  (9) 


and  there  is  no  gas  diffusion  flux  passing  through  the  boundary  of 
fiber  surface: 

Db(B-^j=  0.  (10) 

The  continuity  of  the  diffusive  gas  at  the  boundary  of  the  unit 
cell  results  in: 


and  the  gas  diffusion  rate  is  equal  for  both  at  the  boundary  of  the 
unit  cell: 


Db(B_i) =DeffG-  (12) 

Solving  the  equations  of  (9)— (12),  the  through-plane  effective 
diffusivity  (i.e.  gas  diffusion  normal  to  fiber  arrays  or  layers)  for  both 
the  unit  cell  and  the  whole  ID  regular  fibrous  system  is  obtained: 


Deff  R2~r2 
Db  R2  ;  r2‘ 


(13) 


or  in  forms  of  porosity  c. 


Deff 

Db  ' 


(14) 


Then  we  consider  the  more  complicated  and  realistic  structures,  in 
which  fibers  are  placed  randomly.  Diffusion  process  through  the 
fibrous  system  becomes  more  complex  because  of  the  disorder  of  the 
fiber  arrangement.  In  the  fibrous  system  composed  of  randomly 
located  fibers,  one  fiber  is  assumed  to  be  contained  by  a  polygonal  cell 
whose  boundaries  are  defined  by  the  perpendicular  bisectors  of  the 
lines  joining  each  fiber  with  its  nearest  neighbor,  as  presented  in 
Fig.  2.  The  polygonal  cell  is  called  Voronoi  Tessellation  [28],  which  is 
used  to  characterize  the  randomness  of  fiber  distribution  in  this  study. 

It  is  reasonable  to  assume  that  the  highly  porous  fibrous  system 
is  homogeneous  macroscopically  with  a  constant  concentration 
gradient,  but  different  diffusion  coefficients  in  each  cell.  Thus,  the 
effective  diffusive  flux  Jeff  of  the  system  can  be  obtained  in 
a  volume-averaged  form: 


_  jjdS  _  EJ(M)AS  _  El=0 J(Si)Si 

Jef{  JdS  EM  E»05  ’ 


where  J(S,)  and  S;  are  the  flux  and  area  of  the  ith  unit  cell,  respec¬ 
tively,  N  is  the  number  of  the  unit  cells  in  the  system,  and  (S)  is  the 
mean  area  of  the  unit  cells.  As  N  is  very  large  and  Si  is  very  small 
comparing  with  the  crossing  area  of  the  whole  system,  the  effective 
diffusive  flux  Jeff  can  be  expressed  as  a  function  of  probability 
density  of  unit  cell  area/(S)  in  an  integral  form: 

Jeff  =  J  J(S)/(S)dS.  (16) 

Substituting  Eq.  (16)  into  Eq.  (1 )  leads  to  the  effective  diffusivity: 

Deff  =  J  D(S)f(S)dS.  (17) 


o 
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Fig.  2.  A  unit  cell  in  randomly  distributed  array 
same  area  with  the  circle  in  dotted  line. 


BR  +  -  =  CR, 


(11) 


’  of  parallel  fibers.  The  voronoi  has  the 
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The  local  cell  area  distribution  in  the  Voronoi  Tessellation  is  well 
described  by  the  Gamma  distribution  [29]: 


m 


a^S01-'  _ax 

(S)aWf 


(18) 


For  gas  diffusion  parallel  with  fibers,  the  tortuosity  of  the  ID 
fibrous  medium  is  found  to  be  almost  independent  of  its  porosity, 
and  the  in-plane  effective  diffusivity  (i.e.  gas  diffusion  parallel  with 
fiber  arrays  or  layers)  for  both  regular  and  random  1 D  fibers  is  given 
by  [31]: 


where  F(a)  is  a  Gamma  distribution,  and  a  is  the  scale  parameter.  De ff 

It  satisfies  the  following  equations:  Dj, 


(24) 


J  /(S)dS  =  1,  (19) 

and 

/  S/(S)dS  =  (S).  (20) 


In  analogous  to  the  assumption  successfully  applied  to  the  study 
of  more  structure-sensitive  gas  flow  [24],  this  model  assumes  that 
the  diffusivity  within  the  voronoi  is  identical  to  that  within  the 
circle  of  the  same  area  as  shown  in  Fig.  2.  Substituting  Eq.  (13)  and 
Eq.  (18)  into  Eq.  (17),  the  normalized  through-plane  effective 
diffusivity  of  ID  random  fibrous  media  is  derived: 


Pgff 

Db 


S-nr2  a01  /  S  \“_1  . 
S  +  irr2r(a){(S)J  e 


(21) 


Re-writing  Eq.  (21)  with  t  =  1  -  ( icr2/(S ))  and  x  =  S/(S),  we 
can  obtain  the  following  equation  of  normalized  through-plane 
diffusivity: 


Deff 

Db 


l 


a01  x  -  1  + 

r{a)x  + 1  - 


*-1e-axdx. 


(22) 


The  scale  parameter  a  is  numerically  determined  by  Monte 
Carlo  Simulation.  For  randomly  distributed  fibers,  a  =  3.5  [30]  and 
Eq.  (22)  becomes: 


Deff 

Db 


"  X  —  1  +  £ 


(23) 
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The  model  for  ID  random  fibers  described  by  Equation  (23)  also 
applies  to  2D  random  fibrous  structures,  since  the  through-plane 
effective  diffusivity  is  not  sensitive  to  in-plane  fiber  orientation. 
This  was  demonstrated  from  the  past  numerical  simulation  [14] 
and  can  also  be  inferred  from  the  fact  that  gas  diffusion  is  less 
affected  by  the  presence  of  fibers  than  the  pressure-driven  gas  flow, 
and  the  gas  permeability  was  shown  to  be  independent  of  in-plane 
fiber  orientation  [32,33], 

For  3D  fibrous  materials  with  some  fibers  more  aligned  with  the 
principal  direction  of  gas  diffusion,  the  exact  solution  of  gas  diffu¬ 
sivity  is  difficult  to  obtain.  In  this  study,  the  effective  diffusivity  is 
estimated  from  a  combination  of  local  diffusivities  of  ID  fiber  arrays 
normal  to  and  that  of  ID  fiber  arrays  parallel  with  the  diffusion 
direction  (see  Fig.  3).  Although  the  distribution  of  fibers  in  a  3D 
fibrous  material  is  complex,  it  can  be  considered  that  the  entire  3D 
fibrous  system  is  homogeneous  at  macro-scale  and  the  average 
porosity  of  the  fraction  of  the  fibers  oriented  in  the  three  principle 
directions  (x-,y-  or  z-direction)  shares  the  same  average  porosity  of 
the  whole  fibrous  system.  Three  mixing  rules,  which  were 
proposed  for  permeability  prediction  [25],  are  also  used  in  this 
study  to  determine  effective  diffusivities  based  on  the  mathemat¬ 
ical  analogy  between  Darcy’s  law  and  Fick’s  law.  The  three  mixing 
rules  based  on  different  approximations  and  as  functions  of  fiber 
fraction  q>  =  1  —  c  are  listed  as  follows: 


De  ffW  =  ^DnormW+^DparW, 

(25) 

DeffW  =  [^D^rmW+^Dpiw]  , 

(26) 

DeffW  =  [Dnorm(0)]^[Dpar(0)]^, 

(27) 

||  H  H  U  Diffusion 

Fig.  3.  Illustration  of  3D  or  2D  fibrous  media  composed  of  ID  fiber  arrays. 
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where,  “norm”  and  “par”  mean  normal  to  and  parallel  with  the 
diffusion  direction,  respectively.  In  Eq.  (25),  the  volume-weighted 
diffusivity  assumes  fibers  parallel  to  the  diffusion  direction 
and  fibers  normal  to  the  diffusion  direction  as  parallel  resistors.  In 
Eq.  (26),  fibers  parallel  to  the  diffusion  direction  and  fibers  normal 
to  the  diffusion  direction  are  considered  in  series.  Eq.  (27)  gives 
a  geometric  mean  by  a  mathematical  blend  without  physical 
meaning. 

3.  Model  validation  and  discussions 

Prediction  of  the  proposed  model  is  compared  with  numerical, 
experimental  and  analytical  results  available  in  the  literature.  Since 
the  data  in  the  literature  were  reported  in  terms  of  tortuosity, 
equivalent  diffusivity,  diffusive  resistance,  or  other  related  param¬ 
eters,  they  were  first  converted  to  normalized  effective  diffusivity 
for  direct  comparison  with  the  prediction  of  the  present  model.  For 
the  purpose  of  this  comparison,  we  assume  </>norm  —  0 para  —1/2  for 


in-plane  gas  diffusion  of  2D  fibrous  media;  and  we  assume 
0norm  =  2/3,  and  0para  =1/3  for  through-plane  gas  diffusion  of  3D 
fibrous  media. 

The  normalized  through-plane  effective  diffusivities  calculated 
for  ID  fiber  arrays  in  regular  and  random  packing  are  plotted  in 
Fig.  4.  For  ID  fiber  arrays  in  regular  packing,  the  results  of  the 
present  model  match  almost  perfectly  with  the  numerical  results  of 
Perrins  et  al.  [34]  and  the  experimental  data  from  Kharadly  and 
Jackson  [35],  However,  the  results  of  our  model  are  considerably 
different  from  those  of  Shen  and  Springer’s  model  [17],  who  did  not 
consider  the  width  variation  of  the  void  channel  between  fibers. 
The  prediction  of  Li  et  al.  [19]  is  a  little  higher  than  that  of  the 
present  model,  which  may  arise  from  their  over-idealized 
assumption  that  gas  concentration  is  uniform  at  any  cross-section 
of  the  void  channel. 

Substituting  Eq.  (14)  into  Eq.  (2),  we  can  get  the  tortuosity  of  ID 
fiber  arrays  with  regular  packing  when  the  fibers  are  normal  to  gas 
diffusion: 
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Fig.  8.  Comparison  of  the  mixing  laws  of  2D  in-plane  diffusivity  with  numerical 


(29) 


where  k  is  an  empirical  constant.  In  this  study,  Eq.  (29)  is  employed  to 
determine  the  tortuosity  of  fibrous  media,  and  a  semi-analytical 
model  of  effective  diffusivity  is  obtained  based  on  Eqs.  (2)  and  (29),  viz. 
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Apparently,  Eq.  (30)  can  be  deduced  to  Eq.  (14)  for  ID  regular 
fibers  when  k  =  1.  In  addition,  it  can  be  seen  in  Fig.  4  that  Eq.  (30) 
matches  excellently  with  the  ID  random  model  of  Eq.  (23)  by 
adjusting  the  constant  k  =  1.7.  Therefore,  the  semi-analytical  model 
of  ID  randomly  located  fibers  is  given  by: 
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(31) 


Fig.  4  also  shows  that  the  through-plane  effective  diffusivity  of 
ID  random  structure  is  slightly  lower  than  that  of  ID  regular 
structure,  which  can  be  explained  by  the  fact  that  more  random  and 
disordered  structures  lead  to  longer  tortuous  paths  for  diffusive 
molecules.  When  the  scale  parameter  a  in  Eq.  (22)  increases,  the 
cell  areas  of  fiber  alignments  in  Fig.  2  become  more  uniformly 


distributed.  The  ID  through-plane  diffusivity  from  Eq.  (22)  tends  to 
be  equal  to  that  from  Eq.  (14)  when  a  is  very  large,  as  seen  in  Fig.  5. 
It  is  interesting  to  note  that  the  effect  of  randomness  on  viscous  gas 
flow  though  fibrous  media  is  different  from  the  effect  of  random¬ 
ness  on  gas  diffusion.  The  gas  permeability  of  fibrous  media 
increases  with  increasing  randomness  [37],  but  the  gas  diffusivity 
reduces  with  increasing  randomness  of  fiber  arrangements. 
Therefore,  when  fibers  are  less  randomly  located,  the  fibrous 
structures  will  have  lower  gas  flow  permeability  but  higher  gas 
diffusivity,  making  them  better  candidates  of  breathable  materials. 

In  order  to  validate  the  present  model,  past  numerical  and 
experimental  results  of  GDLs  are  plotted  and  compared  with  the 
prediction  of  the  present  model  in  Fig.  6.  Tomadakis  and  Sotirchos 
[14]  investigated  the  diffusivities  through  ID  and  2D  randomly 
located  fiber  assembles  using  Monte  Carlo  simulation  and  found 
that  ID  and  2D  through-plane  diffusivities  are  similar.  As  can  be 
seen  from  Fig.  6,  both  ID  and  2D  through-plane  diffusivities  from 
Tomadakis  and  Sotirchos  [14]  are  close  to  our  2D  (ID)  through- 
plane  model  of  Eq.  (23),  which  shows  that  ID  and  2D  through- 
plane  diffusivities  are  approximately  independent  of  in-plane 
fiber  orientation.  We  can  also  see  from  the  figure  that  our  model 
prediction  agrees  well  with  the  experimental  data  from  Penman  [4] 
and  Moest  et  al.  [10]  and  the  numerical  simulation  results  by  Nam 
and  Kaviany  [21]  and  by  Becker  et  al.  [22,23],  Nevertheless,  it 
should  be  noted  that  our  model  prediction  deviates  greatly  from 
some  experimental  data  in  the  literature  [9—11],  This  is  due  to  the 
fact  that  chemical  binders  were  applied  in  the  samples  of  these 
experiments,  but  not  considered  in  our  model.  The  binder  used  to 
bind  carbon  fibers  together  can  fill  the  voids  between  the  fibers  (see 
Fig.  7  from  Ref.  [9]).  Therefore,  the  effective  diffusivities  measured 
in  these  experiments  (e.g.  the  experimental  results  of  Toray  TGP- 
060  carbon  fibers  by  Flueckiger  et  al.  [9]  and  Toray  TGP-120 
carbon  fibers  by  Moest  et  al.  [10]  and  LaManna  and  Kandlikar 
[11  ])  are  significantly  lower  than  the  model  prediction.  Moest  et  al. 
[10]  also  measured  effective  diffusivity  of  Freudenberg  FB-A- 
0  carbon  fiber  layers,  which  are  bonded  mechanically  without 
chemical  binder,  and  the  experimental  data  were  comparable  with 
the  present  model  of  2D  through-plane  diffusivity.  In  addition, 
numerical  simulations  by  Becker  et  al.  [38]  show  that  the  virtually 
created  GDL  with  binder  converged  on  the  fiber  contacting  area 
(V0)  had  higher  diffusivity  than  the  reconstructed  GDL  based  on  the 
tomographic  image  of  real  TGP-060  paper  (TO),  as  seen  in  Fig.  6. 
Therefore,  high  performance  of  PEMFCs  with  improved  gas  diffu¬ 
sivity  can  be  achieved  with  less  binder  used  or  less  binder 
distributed  in  void  pores  between  fibers.  We  also  provide  a  semi- 
analytical  model  fitted  by  the  diffusivities  of  TCP  carbon  fibers 
[9-11]  based  on  Eq.  (30): 
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Volume-wersdited  diffusivity  model  pre  g  e  nt  model  Volume-weighted  resistance  model 

Fig.  9.  Relationship  between  the  realistic  model,  volume-averaged  diffusivity  model,  and  volume-averaged  resistance  model. 
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In  Fig.  8,  the  2D  in-plane  diffusivity  of  a  fibrous  material,  con¬ 
sisting  of  ID  random  fibers,  assuming  half  fibers  parallel  and  half 
normal  to  flux  direction,  predicted  by  the  present  model  using 
three  different  mixing  rules,  is  compared  with  the  numerical 
results  from  Tomadakis  and  Sotirchos  [14]  and  Becker  et  al.  [22],  It 
is  interesting  to  see  that  there  is  only  very  small  difference  between 
the  results  based  on  the  three  mixing  rules,  and  the  model 
prediction  using  either  of  the  mixing  rules  agree  well  with  the 
numerical  results.  The  volume-averaged  diffusivity  model  and  the 
volume-averaged  resistance  model  are  the  upper  and  lower  limits 
of  diffusivity  estimates,  respectively  (see  Fig.  9),  while  the 
geometric  mean  is  merely  a  mathematical  estimate  without  phys¬ 
ical  meaning.  Therefore,  a  better  estimate  would  be  the  average  of 
volume-averaged  diffusivity  model  and  volume-averaged  resis¬ 
tance  model,  which  is  adopted  in  this  study  and  plotted  in  Fig.  8: 


—  1 D  through-plane 
•  2D  through-plane 
<  2D  in-plane 
• — 3D  through-plane 
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Based  on  Eqs.  (24),  (31)  and  (33),  the  2D  in-plane  diffusivity 
model  can  be  expressed  as: 
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Likewise,  comparison  of  three  different  mixing  rules  is  carried 
out  for  3D  through-plane  diffusivity  in  Fig.  10.  The  prediction  of  our 
model  for  the  3D  fibrous  structure,  assuming  it  consists  of  ID 
random  fibers  with  1/3  parallel  with  and  2/3  normal  to  diffusion 
direction,  agrees  well  with  the  numerical  data  by  Tomadakis  and 
Sotirchos  [14],  Based  on  Eqs.  (24),  (31)  and  (33),  we  can  derive 
a  semi-analytical  model  for  3D  fibrous  structure: 


Fig.  11.  Effective  diffusivities  for  different  fiber  orientations. 


Fig.  11  compares  the  effective  diffusivities  of  randomly  distrib¬ 
uted  fiber  assembles  with  different  fiber  orientations.  The  ID  in¬ 
plane  diffusivity  is  highest,  which  is  easily  understood  as  its 
tortuosity  is  unity  when  the  fibers  parallel  with  the  diffusion 
direction.  Both  ID  fiber  arrays  and  2D  fibrous  materials  normal  to 
the  diffusion  direction  have  the  highest  diffusive  resistance.  The 
through-plane  diffusivity  of  3D  fibrous  material  is  slightly  lower 
than  the  in-plane  diffusivity  of  2D  fibrous  material,  both  of  which 
fall  between  the  through-plane  and  in-plane  diffusivities  of  ID 
fiber  arrays.  In  general,  the  effective  diffusivity  is  enhanced  with 
increasing  fibers  aligned  with  the  diffusion  direction. 

4.  Conclusions 
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Fig.  10.  Comparison  of  the  mixing  laws  of  3D  through-plane  diffusivity  with  numerical 


Theoretical  models  and  semi-analytical  equations  are  proposed 
to  predict  the  effective  diffusivities  of  ID  fiber  arrays,  and  2D  and 
3D  fibrous  materials  widely  used  as  GDLs  in  PEMFCs.  The  predicted 
effective  diffusivities  agree  well  with  experimentally  measured 
values  and  past  numerical  simulation  results  of  GDLs  containing  no 
chemical  binder.  Specifically,  it  is  found  that  materials  with 
randomly  distributed  fibers  have  lower  diffusivities  than  those  with 
fibers  orderly  distributed.  Using  less  chemical  binders  and 
prevention  of  chemical  binders  filling  the  pores  of  the  fibrous 
materials  increase  effective  diffusivities.  The  present  investigation 
is  significant  to  the  design  and  optimization  of  fibrous  materials  for 
GDLs  and  other  applications. 
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Nomenclature 

a:  edge  length  of  square,  m 

Dh:  bulk  diffusivity,  m2  s-1 
Deft  effective  diffusivity,  m2  s-1 
J:  diffusive  flux,  g  s-1  m-2 
l:  radial  distance  of  a  cell,  m 
r:  fiber  radius,  m 
R:  cell  radius,  m 
S:  cell  area,  m2 

Creek  symbols 

a:  scale  of  gamma  distribution 

fiber  fraction 
c:  porosity 

ep:  percolation  threshold 
8:  azimuth 

r:  gamma  distribution  function 
t:  tortuosity 


